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We model the nonlinear saturation of the r-mode instability via three-mode couplings and the 
effects of the instability on the spin evolution of young neutron stars. We include one mode triplet 
consisting of the r-mode and two near resonant inertial modes that couple to it. We find that the 
spectrum of evolutions is more diverse than previously thought. We start our evolutions with a star 
of temperature ~ 10^" K and a spin frequency close to the Kepler break-up frequency. We assume 
that hyperon bulk viscosity dominates at high temperatures (T ~ 10^ — 10^° K) and boundary 
layer viscosity dominates at lower temperatures (~ a few x 10* K). To explore possible nonlinear 
behavior, we vary properties of the star such as the hyperon superfluid transition temperature, the 
strength of the boundary layer viscosity, and the fraction of the star that cools via direct URCA 
reactions. The evolution of the star is dynamic and initially dominated by fast neutrino cooling. 
Nonlinear effects become important when the r-mode amplitude grows above its first parametric 
instability threshold. The balance between neutrino cooling and viscous heating plays an important 
role in the evolution. Depending on the initial r-mode amplitude, and on the strength of the viscosity 
and of the cooling this balance can occur at different temperatures. If thermal equilibrium occurs 
on the r-mode stability curve, where gravitational driving equals viscous damping, the evolution 
may be adequately described by a one-mode model. Otherwise, nonlinear effects are important 
and lead to various more complicated scenarios. Once thermal balance occurs, the star spins- 
down oscillating between thermal equilibrium states until the instability is no longer active. The 
average evolution of the mode amplitudes can be approximated by quasi-stationary states that are 
determined algebraically. For lower viscosity we observe runaway behavior in which the r-mode 
amplitude passes several parametric instability thresholds. In this case more modes need to be 
included to model the evolution accurately. In the most optimistic case, we find that gravitational 
radiation from the r-mode instability in a very young, fast spinning neutron star within about 1 
Mpc of Earth may be detectable by advanced LIGO for years, and perhaps decades, after formation. 
Details regarding the amplitude and duration of the emission depend on the internal dissipation of 
the modes of the star, which would be probed by such detections. 

PACS numbers: 04.40.Dg, 04.30.Db, 97.10.Sj, 97.60.Jd 



I. INTRODUCTION 

Neutron stars are believed to be born in the aftermath 
of core-collapse supernova explosions as the stellar rem- 
nant becomes gravitationally decoupled from the stellar 
ejecta. Two interesting and timely question are: (1) Do 
neutron stars spin at millisecond periods at birth or do 
they spin closer to the observed periods of young pul- 
sars?; (2) Do they emit gravitational radiation that is 
detectable by interferometers on Earth? 

Theoretically, conservation of angular momentum in 
the core collapse of a 8 — 30Mq progenitor can lead to a 
newborn neutron star with a period of ~ I ms or shorter. 
Observationally, the fastest known young pulsar is in the 
Large Magellanic Cloud supernova remnant N157B, and 
has a rotation period of 16 ms f]\. The Crab pulsar is the 
next fastest neutron star in a supernova remnant with an 
age ~ 10^ yr. Its current rotation period is 33 ms. Its 
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initial period is estimated to be ~ 19 ms j2j by assuming 
the rotational spin-down is well described by a power law 
fl cx -17" with braking index n = 2.51 ± O.OI Q. One 
way to predict the distribution of initial pulsar periods is 
through population synthesis studies. These studies gen- 
erally use present day observations with some assumption 
of their time evolution to reconstruct the birth distribu- 
tion periods and magnetic fields of the pulsar population. 
Current studies favor initial periods in the range of sev- 
eral tens to several hundreds of milliseconds [4, 5] . The 
apparent discrepancy between the theoretically possible 
fast rotation rates and the observed slow rotation rates 
of young neutron stars could be reconciled if the r-mode 
instability or some other mechanism could spin neutron 
stars down efficiently, preventing them from maintaining 
millisecond periods. 

R-modes are quasi-toroidal oscillations in rotating flu- 
ids that occur because of the Coriolis effect. These modes 
are driven unstable by gravitational radiation reaction 
via the Chandrasekhar-Friedman-Schutz (CFS) mecha- 
nism d, 0]- In the absence of fluid dissipation, the CFS 
mechanism causes any mode that is retrograde in the 
co-rotating frame, but prograde in the inertial frame, to 
grow as gravitational radiation is emitted [1, 0]. The 
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most unstable r-mode is the n = 3, m = 2 mode, called 
the r-mode throughout the rest of the paper (Here n and 
m label the degree and order of the Legendre functions 
associated with the mode.) The gravitational driving 
equals viscous dissipation for this mode along a critical 
curve in the angular velocity - temperature (f2 — T) phase 
space. Above the critical curve the r-mode is linearly 
unstable and the r-mode amplitude grows exponentially. 
Once the r-mode amplitude passes its first parametric 
instability threshold, other near-resonant inertial modes 
are excited via energy transfer from the r-mode. At this 
point nonlinear effects become important. The paramet- 
ric instability threshold amplitude depends on internal 
neutron star physics and changes with angular velocity 
and temperature. 

The first investigation that modeled the spin-down of a 
young neutron star due to the r-mode instability was per- 
formed by Lindblom et al. fT3| (See also Owen et al. [ll| 
for a more detailed analysis.) They used a simple one- 
mode evolution model that assumes the r-mode ampli- 
tude saturates because of nonlinear effects at some arbi- 
trarily fixed value. The saturation amplitude was chosen 
to be of order 1. In this model, once the instability is sat- 
urated, the star spins down at fixed r-mode amplitude. 
Lindblom et al. estimated that a newborn neutron star 
would cool to approximately 10^ K and spin down from a 
frequency close to the Kepler frequency to about 100 Hz 
in ^ 1 yr. In their calculation they included the effects 
of shear viscosity and bulk viscosity for ordinary neutron 
star matter composed of neutrons, protons and electrons 
and assumed modified URCA coohng. 

Jones [13] and Lindblom & Owen [T^ pointed out that 
if the star contains exotic particles such as hyperons, in- 
ternal processes could lead to a very high bulk viscosity 
in the cores of neutron stars. They predicted that for 
young neutron stars this viscosity would either eliminate 
the instability altogether or leave a short window of in- 
stability of up to a day or so for modified URCA cooling 
p^ . that would render the gravitational radiation un- 
detectable. Andersson, Jones, and Kokkotas [l8| found 
that, in the case of strange stars, young neutron stars can 
evolve along the r-mode stability curve reaching a quasi- 
adiabatic equilibrium at low r-mode amplitudes. Similar 
conclusions were reached by Reisenegger & Bonacic 
for large hyperon bulk viscosity. We include nonlinear 
effects and find that the spectrum of possible evolutions 
is more diverse. Quasi-equilibrium evolutions along the 
r-mode stability curve are just one scenario. 

Schenk et al. [lO*! developed a formalism to study 
nonlinear interactions of the r-mode with other inertial 
modes. They assumed a small r-mode amplitude and 
treated the oscillations of the modes with weakly nonlin- 
ear perturbation theory via three-mode couplings. This 
assum ptio n was tested by Arras et al. [20] and Brink 
et al. |2ll . [2^ . m, [23] . Arras et al. proposed that a tur- 
bulent cascade will develop in the strong driving regime. 
They estimated that the r-mode amplitude was small and 
could have values between ID"**— 10~^. Brink et al. com- 



puted the interactions of about 5000 modes via approx- 
imately 1.3 X 10^ couplings among modes with n < 30. 
They modeled the star as incompressible and calculated 
the coupling coefficients analytically. The couplings were 
restricted to near resonant modes with a fractional de- 
tuning of Su!/{2n) < 0.002. Brink et al. showed that the 
nonlinear evolution saturates at amplitudes comparable 
with the lowest parametric instability threshold. They 
did not include spin or temperature evolution in their 
calculation. 

In a previous paper [25j we investigated the nonlin- 
ear saturation of the r-mode instability for neutron stars 
in LMXBs and included temperature and spin evolution. 
In this paper we begin a study of the nonlinear develop- 
ment of the r-mode instability in newborn neutron stars. 
We use an effective three-mode treatment in that we in- 
clude one triplet of modes with a statistically relevant 
coupling and detuning coefficient and treat it as the low- 
est parametric instability threshold. The triplet consists 
of the n = 3, m = 2 r-mode and two near-resonant in- 
ertial modes with n = 14 and n = 15 that couple to 
it. The exact inertial modes that are excited when the 
r-mode grows above its parametric instability threshold 
will change as the star spins down. However, since the 
nearest resonance is statistically expected to be at a de- 
tuning of 5uj/2n 10~^ and the mode coupling coef- 
ficient for the lowest parametric instability threshold is 
typically of order one |2j], a mode triplet using these 
values should provide a qualitatively correct picture. 

The model also includes neutrino cooling via a combi- 
nation of fast and slow processes, viscous heating due 
to hyperon bulk viscosity and boundary layer viscos- 
ity, and spin-down due to gravitational radiation and 
magnetic dipole radiation. In oder to explore possible 
nonlinear behaviors we vary: (1) the hyperon super- 
fluid transition temperature Th, which is believed to be 

10^ — 10^° K; (2) the strength of the hyperon bulk vis- 
cosity; (3) the boundary layer viscosity via the slippage 
factor Sns that parametrizes the interaction between os- 
cillating fluid core and the elastic crust [H, [13] ; and (4) 
the fraction of the star that cools via direct URCA reac- 
tions /dU- 

We find a variety of scenarios that depend on these 
parameters and on the initial r-mode amplitude. The 
star cools until the cooling is balanced by viscous heat- 
ing. It then follows a quasi-adiabatic evolution either on 
the r-mode stability curve or on other Cooling = Heat- 
ing curves on which the neutrino cooling is balanced by 
viscous heating from the three modes. For low hyperon 
bulk viscosity or when we include only slow cooling, we 
find runaway behavior in which the energy dissipated by 
the two inertial modes is not sufficient to stop the growth 
of the r-mode amplitude and several parametric instabil- 
ity thresholds are passed. Modeling such behavior ac- 
curately requires the inclusion of multiple mode triplets, 
which is beyond the scope of this paper. 

The three mode problem is a natural first step in study- 
ing nonlinear effects for the r-mode instability. This 



model is valid as long as the r-mode amplitude does 
not grow above parametric instability thresholds that 
are higher than the lowest threshold. The thresholds 
are functions of the angular velocity and temperature 
of the star and depend on the viscous damping rates of 
the modes and model details. The star cools at constant 
angular velocity until the cooling is stopped by viscous 
heating and then spins down oscillating between thermal 
equilibrium states until the r-mode is no longer unsta- 
ble. The r-mode amplitudes we find are still fairly low 
^ 10^^ — 10^'^ and the spin-down torque due to grav- 
itational radiation reaction is typically lower than that 
due to magnetic dipole radiation. So, the spin down 
timescale, which is approximately the timescale on which 
the r-mode is unstable, is dominated by the magnetic 
dipole radiation timescale for B ~ 10^'^ G. Our evolu- 
tions are determined by competitions between cooling 
and heating, gravitational driving and viscous damping, 
and magnetic dipole and gravitational spin-down with 
the competition between cooling and heating playing the 
most prominent role. We expect that some of this be- 
havior will hold for more sophisticated models as well. 

Our evolutions are determined by competitions be- 
tween neutrino cooling and viscous heating, gravita- 
tional driving and viscous dissipation, neutrino cooling 
and gravitational driving, neutrino cooling and magnetic 
spin-down of the star with the competition between cool- 
ing and heating placing the most prominent role. The 
three-mode model is adequate as long as the mode am- 
plitudes do not grow high enough to excite more modes 
in the star. Stable three-mode evolutions are more likely 
to happen for high viscosity such that due hyperons. 

The results are summarized in Sec. |lTl Model details 
are presented in Sec. IIIIl Stable evolutions are examined 
in detail Sec. IIVI The runaway evolutions are discussed 
in Sec. |Vl The prospects of detecting gravitational radi- 
ation are considered in Sec. lVIl Limitations of the model 
are discussed in Sec. I VIII Concluding remarks are pre- 
sented in Sec. IVlIIl 



II. SUMMARY OF RESULTS 

Our evolutions start with an initial temperature of 10^*^ 
K and an initial rotation rate at the Kepler break-up fre- 
quency fl = 0.67-y/7rGp, where p is the average density of 
the star. These initial conditions place the neutron star 
above the r-mode stability curve in the unstable region. 
The evolutions are independent of the initial tempera- 
ture for temperatures for which the cooling rate is faster 
than the gravitational driving rate. The star cools at ap- 
proximately constant angular velocity until the cooling is 
stopped by viscous heating. Once thermal equilibrium is 
reached, the star spins down oscillating though a series of 
thermal equilibrium states. The spin-down is the longest 
part of the evolution lasting at least several years, while 
the initial cooling at constant angular velocity lasts min- 
utes to hours depending on whether low or fast cooling 
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FIG. 1: Schematic plot showing different trajectories that 
the star can follow in the Cl - Tg plane. Here Cl — 
Q.C = 0FG^, and Tg = T/{10^K). The Cooling = Heating 
(C = H) curves begin when the star starts spinning down and 
end when it reaches the stability curve. For simplicity we did 
not include mixed type I-H and type II-IH evolutions, where 
the thermal equilibrium becomes unstable and the star cools 
from one C = H curve to another. A flow chart that includes 
all possible scenarios is displayed in Fig. |2l 



processes dominate. So, the timescale on which the r- 
mode instability is active and gravitational radiation is 
emitted is approximately the spin-down timescale. 

The spin-down part of the trajectory of the star can be 
approximated by a Cooling — Heating (C — H) curve. 
This curve is not unique and is determined by equat- 
ing the neutrino cooling and the viscous heating using 
quasi-stationary amplitude solutions (Eqs. or Eq. 
(dni). Our schematic hyperon bulk viscosity model cor- 
responds to an r-mode stability curve with one maximum 
(See Fig. [1]). We find three possible C — H curves: one 
on each side of the maximum, and one that is on the r- 
mode stability curve, where gravitational driving equals 
viscous dissipation. The latter corresponds to evolutions 
for which a one-mode model is adequate and the former 
two C = H curves approximate triplet evolutions. 

The trajectory that the star follows depends on 
whether the r-mode amplitude passes its parametric in- 
stability threshold before the r-mode stability curve or 
afterwards. If the parametric instability threshold is 
reached first, the star will spin down on the C — H near 
the high T side (T > Tpoak) of the r-mode stability curve 
(type I evolution). In this scenario, the viscous heat- 
ing due to all three modes balances the cooling. If the 
r-mode stability curve is reached first, the star can spin 
down following the r-mode stability curve (type II evolu- 
tion) or, for lower initial r-mode amplitude, cool through 
that thermal equilibrium region and spin down on the 
C = H near the low T side (T < Tpeak) of the r-mode 
stability curve (type III evolution). As the star spins 
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FIG. 2: Flow diagram for the various evolution scenarios. The parametric instability threshold is abbreviated as PIT. 



down and cools slowly, the thermal equilibrium can be- 
come unstable. When this happens the star cools to the 
next closest thermal equilibrium region leading to mixed 
type I-II (See Fig. [5] for an example.) and type II-III 
(See Fig. [8]) evolutions. 

The quasi-stationary solutions are independent of ini- 
tial conditions. Evolutions in which the three modes are 
sufficient to model the behavior follow trajectories that 
oscillate around C — H curves. So, they are roughly in- 
dependent of initial conditions as well. However, the ini- 
tial r-mode amplitude can determine which type of evo- 
lution will occur. At a given viscosity, cooling and initial 
spin- frequency, for initial r-mode amplitudes above a cer- 
tain Ca (0) the r-mode amplitude will always reach para- 



metric instability before the stability curve and hence a 
type I evolution will occur. And below this value type 
II evolutions will occur. A similar amplitude threshold 
exists between type II and type III evolutions. 

If no thermal equilibrium can be reached, the r-mode 
amplitude grows above several parametric instability 
thresholds and excites more inertial modes. In order to 
accurately model runaway behavior more mode triplets 
need to be included. This is beyond the scope of this 
paper. We obtain runaway scenarios by either lowering 
the hyperon bulk viscosity or by turning off the fast neu- 
trino cooling. So, for low mass neutron stars that have 
lower central densities and no hyperon population, multi- 
ple triplets of modes are required to model the nonlinear 
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behavior accurately. If energy transfer to inertial modes 
does not stop the instabihty at a low r-mode amplitude, 
other nonlinear saturation effects would become impor- 
tant such as nonlinear bulk viscosity [l5|. 

The evolutionary scenarios described here are summa- 
rized in Fig. [21 Each stable evolution scenario is pre- 
sented in more detail in Sec. IIVI and the unstable evo- 
lutions are discussed in Sec. El We consider an evolu- 
tion to be 'stable' when one mode triplet is adequate 
to model the evolution. In this case the neutrino cool- 
ing is stopped by viscous heating due to the three modes 
and the star spins down oscillating between thermal equi- 
librium states. 'Unstable' or 'runaway' evolutions occur 
when the star reaches thermal equilibrium only after the 
r-mode and the inertial modes overshoot several para- 
metric instability thresholds reaching amplitudes of or- 
der unity. Modeling such scenarios accurately requires a 
larger mode network. 



III. SETUP 



A. Three-mode Evolution Equations 



AofRef. [2|. 

We further rescale the equations: 
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dC-y 
~d^ 



|ow| \5uj\Q. 



If} 



\6oj\n 



2V^ 

i 



y=Cf3Cj, (2) 



lUJ. 



—J-C 



77 



i 



where we have used the same rescaling as in Ref. [25| with 



7j/r2c, and Sic is a fixed angular 



frequency chosen for reference. 

The rotational phase t is rescaled by the fractional 
detuning as f — T|(5a)| and the mode amplitudes scaled 
by the zero-viscosity parametric instability threshold 
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|(5a)|\/f^ 1^ I |(5<D|%/n^ 

—-j==-, y^m = 

|(5j)|^/n^ 



(3) 



The evolution equations for the mode amplitudes when 
the n = 3, m = 2 r-mode couples to two near-resonant 
inertial modes can be written as 



dC^ 

dCn 
dT 
dC^ 
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The scaled frequency u)j is u)j = tOj/fl, the dissipation 
rates of the inertial modes, also called daughter modes, 
are 7^ and 7^ , 7q is the sum of the driving and damping 
rates of the r-mode 7^ — ^gr ^ lav, and the dimension- 
less coupling IS K, = n/ (MR^ft^). The rotational phase r 
is defined hy dr — dt. In terms of the amplitude vari- 
ables of Schenk et al. [l| and Brink et al. [HI, [H HI 
Cj — ^JVl{t)cj{t), which are normalized to unit energy, 
i.e., mode amplitudes of Cj — 1 correspond to a mode 
energy Ej ^ e ^ MR'^fl'^ equal to the rotational energy 
of the star^. The equations of motion in Schenk et al. 
assume constant S7. We are interested in a situation in 
which the uniform angular velocity of the star changes 
very slowly on a timescale of a rotation period. Writing 
the equation of motion in terms of rotational phase and 
the amplitude variables Cj removes the time dependence 
(Without rescaling, the amplitude evolution equations 
would have terms proportional to dQ/dt.) A derivation 
of the equations of motion for the three mode system 
from the Lagrangian density can be found in Appendix 



B. Driving and Damping Rates 



For our benchmark calculations, we adopt the neutron 
star model of Owen et al. Rcf. [ll| (n = 1 polytrope, 
M = 1.4M0, flc = 8.4 X 10^ rad sec'^ and R = 12.5 km) 
and use their gravitational driving rate for the r-mode 



O.Obsec-^MiAR 
3.26 



4 

12.5 



IkHz 
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sec 



where Mi. 4 = M/I.4M0, i?i2.5 = i?/12.5km and = 
n/ric. If hyperons are present in the core of the neutron 
star then hyperon bulk viscosity dominates other forms 
of bulk viscosity. In this paper we consider hyperon bulk 
viscosity and boundary layer viscosity to be the main 
sources of dissipation. We use a schematic model with a 
few adjustable parameters. See Appendix B for a detailed 
description. 



[2] The conversion between the ampHtudes of Schenk et al. and the 
a used by Lindblom et ah [lo| is Cq = « Ri O.lcj for an 

n = 1 polytrope with J = 1/(MR'^) drpr^. 
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Angular Momentum and Temperature 
Evolution 



Angular momentum is lost via gravitational wave emis- 
sion 



dJ ^ ^ in 

-77 — ^IGbJc inwdc : 

at TM 



(5) 



where Jc rmodo = -ima/uja)ea\ca\'^ ^ -3Mi?2r2|cQ.p = 
— 3Afi?^|CQp, / is the moment of inertia of the star, 
m/TM is the magnetic breaking torque, and tm is the 
corresponding timescale. We have adopted the simplest 
magnetic dipole model with a timescale 
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3 X 109 sec MiaR^ 



(6) 
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6.7 X 108 



where B12 = 5/(10^^ G) and the magnetic dipole mo- 
ment = 10'^''/i3o G cm'^. Fallback disks observed around 
isolated highly magnetic neutron stars [1^ and gas inside 
the magnetosphere could make the spin-down more com- 
plicated. 

Eq. ([5]) can be rewritten in terms of the scaled variables 
in Eq. © as 



dJ 



67GR MR^n,\5Cj\ 



{AkYujisCj'y 
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\5uj\tm 



We then consider the thermal evolution of the star 



(7) 



(8) 



We use a specific heat C{T) = 1.5 x lO^^Tg erg K'^ 
which assumes normal neutrons. Yakovlev and Pethick 
(4l| found that cooling observations are consistent with 
weak neutron superfluidity with critical temperatures of 
r„ < 2 x 10^ K and strong proton superfluidity. The 
energy gap calculations for neutron superfluidity are still 
very uncertain. Our temperatures are typically above 
5 X 10* K and we choose to ignore neutron superfluidity. 
If neutrons are superfluid, the heat capacity is reduced. A 
lower bound to C(T) is the electron contribution, which 
is about 20 times lower than the neutron-neutron specific 
heat used here [i^ . 

The neutrino luminosity is a superposition of fast and 
slow processes 

Lv = L^\jTqR^\j{T /Tp) + Li-avT^Ri-aviT /Tp) (9) 
where the proton superfluid reduction factors for the 



modified and direct URCA reactions are taken from Eqs. 
(32) and (51) of Ref. [13] ■ The other constants are 
adopted from Ref. [13: imU = 1-0 x 10^^° erg sec" \ 
^du = /du X lO'*^ ergsec"^, and the electron-ion and 
neutron-neutron neutrino bremsstrahlung are given by 
Le-i = 9.1 X 10^^ ergsec"^, L^-n ~ 10^* erg sec" ^ 
[3^ . |40| . [4^ . The fraction of the star /du that is above the 
density threshold for direct URCA reactions is in general 
dependent on the equation of state [4l|. We take it to 
be 10% for all stable evolutions in this paper. We use a 
constant Tp — 5.0 x 10® K, which is consistent with Fig. 5 
of A more realistic model would involve a critical 
temperature that changes with the density of the star. 
In terms of the scaled variables Eq. ([5]) becomes 

C(T)— = ; ' (c^a7a.|CaP+^/37/3|C/3p(10) 



■a;-y7^|C^n 



D. Quasi- Adiabatic Spin and Temperature 
Evolution 



In type I and type III evolutions, after an initial pre- 
cursor, the mode amplitudes settle in and later oscillate 
around their quasi-stationary solutions 
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7/3 + 77 - 7a 
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Appendix A provides a brief derivation. 

In the limit of \Suj\/jj << 1, the quasi-steady mode 
amplitudes simplify to 



47a 77 



In the limit when |(5aj|/7j >> 1 

47/377^ l/S |2 ^ 



47a 7/3 

n\scj\^ 



(12) 



I a 



I a 



(7/3+77-7a)^^ 
|2 _ 47a7/3» 

(7/3 + 77 - 7a)^ ' 



47^7^0 



(7/3+77-7")^ 



(13) 



Assuming that dCi/df k, J k, IVt and using 
Eqs. pT|) in Eqs. ([71 fTU|) the spin and thermal evolution 
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equations can be rewritten as 



_ 67GR 7/377 
df (l'^\Sw\ Ak^iiupLu^ 

where / = //(Mi?2), and 



tan^ 



TMflc\6u}\ 

(14) 



7a 7/377 



1 



7a 



+ + 



tan"^ 



(15) 



By setting the right hand side of the above equation to 
zero, one can find the CooUng = Heating (C = H) curve. 
Below, we find that Eqs. (|14 p - (fT5|) approximate well the 
spin-down part of the evolution for type I and III scenar- 
ios throughout the stable regime. 

In type II evolutions the trajectory settles on or os- 
cillates around the low T branch of the r-mode stability 
curve. The thermal equilibrium on the stability curve 
corresponds to zero daughter mode amplitudes and an 
r-mode amplitude determined by the one-mode quasi- 
stationary state 



IG 



C— H one — model 
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2MR'^niVL\5Cj\^ 

(16) 



E. Validity of Approximations 

1. Density of Resonances and the Effective Three Mode 
Coupling Approximation 

We choose a mode triplet with detuning, coupling co- 
efficient, and viscous damping rates that were found to 
be representative for the first near-resonance of a dense 
set of modes confined to the \—2Vt,2Vl] frequency range. 
Brink (See Chapter 5.2 of Ref. [2J|) estimated statisti- 
cally the expected lowest detuning for modes that are 
uniformly distributed. She considered all modes that 
couple to the r-mode and found a detuning of 5u}/{2Vl) ^ 
10^** for n < 14. For a uniform frequency distribution 
this lowest detuning is inversely proportional to the num- 
ber of direct couplings to the r-mode N, which increases 
drastically with n {N « '^max/S for n < ri^ax)- How- 
ever, the viscosity also increases with n. Hyperon bulk 
viscosity dominates at T ~ 10^ K, and we estimate that 
it scales roughly as 7hb oc n^. A very rough upper bound 
of n 100 can be obtained for 7hb ~ tj. The exact n 
cutoff beyond which higher order modes can be ignored 
is model dependent, and finding it involves an extensive 
calculation of damping rates, which is left for future in- 
vestigation. 

Preliminary calculations [43] show that when second 
order corrections to the mode frequencies are included, 
the average lowest detuning is still of order 6u!/{2n) ~ 




FIG. 3: The trajectories in the f2 — Tg plane for evolutions 
with diflterent coupling coefficients are shown for three values 
of k: 3.95, 1.25 and 0.395. All other parameters are the same. 
The hyperon superfluidity temperature is Th — 6.0 x lO'' K. 
The cooling is stopped by the heating later for larger k and 
the thermal oscillations are smaller. 



10^"*. The exact n and m values of the modes that have 
this low detuning change with angular velocity. So, in 
principle, a more accurate modeling of the spin-down pro- 
cess would include multiple mode triplets and allow the 
identities of the modes comprising the lowest parametric 
instability threshold to change as the star evolves. How- 
ever, since the lowest resonance always occurs for similar 
high n modes with damping rates varying within a fac- 
tor of 2, we conjecture that one mode triplet with these 
typical coupling k ~ 1 and detuning coefficients can act 
as an 'effective' lowest resonance. 

To test the dependence of our evolutions on the cou- 
pling coefficient we vary it by factors of -x/TO between 
K = 0.395 and k = 3.95, while keeping all other param- 
eters the same. Fig. [3] shows that the evolutions change 
but are qualitatively similar. Increasing k decreases the 
parametric instability threshold and the quasi-stationary 
mode amplitudes. This lowers the viscous heating due 
to the three modes. The cooling is stopped later increas- 
ing the angular velocity at which the star intersects the 
r-mode stability curve. The system is less dynamic for 
larger k and the amplitudes of the thermal oscillations 
are lower. We also varied the detuning by a factor of 
100, between 10~^ and 10~^, and did not observe sig- 
nificant changes in the evolution. The exact tempera- 
ture at which the star spins down and the initial r-mode 
amplitude |Cq|(0) at which the transition between dif- 
ferent types of evolutions occurs change with k and Suj, 
but, overall, the evolution scenarios we obtain are quali- 
tatively robust. 

A shortcoming of the three-mode model is that there 
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are cases when the r-mode grows above several paramet- 
ric instability thresholds and the system runs away. In 
order to model this scenario accurately one would have to 
include an oscillator network that takes into account not 
only direct couplings to the r-mode, but also couplings 
among inertial modes similar to the work of Brink et al. 




2. Buoyancy Effects 

Buoyancy effects are expected to be important when 
the angular velocity of the star is comparable with the 
Brunt- Vaisala frequency. In a recent paper, Passamonti 
et al. [m compared the inertial modes of a non-stratified 
star with the g-modes of a stratified system. They found 
that each g-mode they studied approached a particular 
inertial mode as the rotation rate was increased. The 
n = 3, m = 2 r-mode does not have a g-mode coun- 
terpart. The conjecture is that if composition gradients 
are included, inertial modes will approach g-modes at 
low rotation rates, when buoyancy dominates. As Q de- 
creases the inertial mode frequencies level off rather than 
decreasing oc fl. Once this happens, the r-mode loses 
its resonant modes, and at low rotation rates the r-mode 
amplitude would grow unabated as long as gravitational 
driving dominates viscous dissipation. 

The Brunt- Vaisala angular frequency is 

where the proton fraction is taken to be cc ~ 5.6 x 
10"'^(p/pnuc)and the nuclear density pnuc = 2.8 x 
10^** gcm~'^ [36|. This corresponds to a spin frequency of 

i^N ~ 160 Hz (18) 

When 2il is of the order N or lower the spectrum of 
inertial modes will change and the r-mode will lose it 
near-resonant modes. However, at « N/2 or « 80 
Hz, the r-mode is in the stable regime for Sns > 0.01 
for our boundary layer viscosity model, and hence this 
window of no nonlinear effects could be irrelevant. 

For our viscosity model for a critical temperature of 
Th = 1.2 X 10^" K, the minimum spin frequency on the 
r-mode stability curve for T < Tpoak and 5ns ~ 0.3 is at 
i^min ~ 300 Hz (51 « 0.22^TrGp). Passamonti et al. find 
that Coriolis effects dominate beyond f] 0.3 [355. 
ignoring buoyancy in our evolutions should not lead to 
large qualitative changes. This frequency scales as 

/ q \ 4/11 / T \ "2/" 
.„,„. 300H.(^) [jj;^) , («) 

when the temperature at the minimum is much smaller 
than the critical temperature for hyperons T^i^ << T^. 
It is important to note that fmin is not very sensitive to 



the critical temperature Th, but moderately sensitive to 
the shppage factor S^^- At = 1.2 x IQi" K we find: 

l^min ~ 90 Hz for S'ns = 0.01, Vm\n ~ 70 Hz for Sns = 

0.005, in rough agreement with the estimated scaling. 
For Th < 10" K and S'ns > 0.05 we obtain v^^in > 50 Hz. 



3. Magnetic Fteld Effects 

The frequency at which the magnetic dipole and grav- 
itational radiation spin-down rates are equal is 

i'kHzlCap/^ « Q.QMI^l'^Bli^. (20) 

Higher ^'kHzlCcp/'^ implies spindown is mainly via gravi- 
tational radiation rather than magnetic dipole radiation, 
and vice-versa. In terms of the scaled r-mode amplitude 
the equation becomes 

« l.\Bli''M~t'^Rl^,^ X (21) 

/ K N,2/3 / 2X 10-n '^' / Cj,3[b^ \ 

VI. 25/ V ) V0-38 x 0.28j ' 

where Ca{^) is a function of fl with the exact angular 
velocity dependence changing with the type of evolution. 
The magnetic dipole spin-down torque becomes larger 
than the gravitational driving torque for IC^I < 4 at 
V < 600Hz and B — 10^^ G. While our simulations reach 
r-mode amplitudes as high as \Ca\ ~ 20, the amount of 
time spent at high r-mode amplitudes is low and at B = 
10^"^ G the spin-down timescale is dominated by dipole 
spin-down for most evolutions. An r-mode amplitude of 
ICqI < 4 corresponds to Cq < 4 x 10"^ and a Schenk et 
al [19] amplitude of |c„| = \Ca\/V^ < 7 x lO"'^. 

Typically, for B ~ 10^^ G we obtain r-mode instability 
timescales of a few years. For B ^ 10^^ G the spin-down 
timescale grows to ~ 10 years or more. This timescale is 
also dependent on the frequency at which the spin-down 
stops, which is partially determined by the shape of the 
r-mode stability curve. 

For low magnetic fields B < 10^^ G the gravitational 
radiation torque is in general larger than the magnetic 
dipole torque and the spin-down timescale ranges from 
tens to thousands of years depending on the strength 
of the viscosity and of the cooling. Since there are no 
observed young millisecond pulsars, scenarios with very 
long spin-down timescales are probably ruled out by ob- 
servations. 

The splitting of the frequencies because of the Lorentz 
force introduced by the magnetic field is small. The cor- 
rections for the inertial modes are about three orders of 
magnitude larger than for the r-mode, but are smaller 
than the typical detuning of 10~^ for moderate magnetic 
fields of B - 10^2 - 1013 G and n < 15 (See Appendix 
C.) We include dipole spin-down, but ignore these small 
frequency corrections. 

We also do not include differential rotation, which, if 



present, could change the modes of the star significantly 

IV. STABLE EVOLUTIONS 

In this section we examine different types of stable 
evolutions in more detail. We assume a slippage fac- 
tor 5„s ~ 0.3, a fraction of the star above the density 
threshold for direct URCA reactions /du = 0.1 and a 
critical temperature in the stellar core for proton su- 
perfluidity Tp = 5.0 X 10^ K. The critical temperature 
for hyperon superfluidity Th is treated as a free param- 
eter and varied between 10^ and 10^° K. We take a 
coupling coefficients of k = 1.25, a fractional detuning 
\6Cj\ = |(5a;|/ri = 2.0 x 10^"*, and the daughter mode 
frequencies tD^ = 0.38, ui^ — 0.28 and viscous damping 
corresponding to the n = 14, m — —5 and n — 15, m — 3 
modes. 



A. Type I 

In this subsection we present two evolutions with hy- 
peron superfluidity temperatures of Th = 2.0 x 10^ K 
(Fig.gD and = 6.0 x 10^ K (Fig.E]). The initial angu- 
lar velocity fli — 0.67 and temperature Tgj = 10 are the 
same in both simulations. The solid arrows in the dia- 
gram in Fig. [2] map the stages that the star goes through 
for the two trajectories. 

The first evolution with = 2.0 x 10^ K is of type I: 
the r-mode amplitude passes its parametric instability 
threshold before reaching the reaching the r-mode stabil- 
ity curve and the star spins down on the C = H curve 
near the high T side of the r-mode stability curve. Fig. 
mja) shows the trajectory in the 51 — T plane. The evolu- 
tion starts in the unstable regime where ^gr » Jav- 
The r-mode amplitude increases exponentially until it 
passes the first parametric instability threshold and ex- 
cites two inertial modes. After a brief precursor of 
about 15 min the three modes settle in their quasi- 
stationary states. The star continues to cool at constant 
angular velocity until it reaches thermal equilibrium at 
Tg « 1.12. This occurs far from the r-mode stability 
curve T9CFs(f^ = 0.67) = 0.47. The neutron star then 
spins down and cools slowly oscillating around the C = H 
curve determined by the quasi-stationary mode ampli- 
tudes. Fig. [31[b) shows the r-mode amplitude oscillating 
around its quasi-stationary solution. The star intersects 
the r-mode stability curve at Jib = 0.27 and Tgb = 0.76. 
The spin-down timescale ia^b ~ 1-1 yrs aX, B = 10^^ G. 
This timescale is independent of detuning, but changes 
with B. For B = 10" G, ta^b « 14.3 yr. Lowering the 
magnetic field further leaves the timescale unchanged be- 
cause beyond this B the gravitational radiation torque is 
larger than the magnetic dipole torque. Point 6 can be 
determined analytically by finding the intersection be- 
tween the r-mode stability curve and the C = H curve. 
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FIG. 4: Type I evolution. Th = 2.0 x lO'-* K. (a) The 
trajectory of the star is shown in the Cl — Tg plane. The star 
cools from Tgi = 10, f2i = 0.67 to Tga ^ 1.12 in about 10 
min. At this point the viscous heating is large enough to stop 
the neutrino cooling. The star then spins down and cools 
to fib ^ 0.27, Tgb « 0.76 in ta^b ~ 1-14 yr. At point b 
the star intersects the r-mode stability curve and the r-mode 
becomes stable. The star starts cooling at constant angular 
velocity until the cooling rate decreases below the magnetic 
dipole spin-down rate, (b) The r-mode and its parametric 
instability threshold amplitude are shown as a function of 
time. As before it can be seen that quasi (Eq. I12|l is a 
good approximation for the average r-mode amplitude until 
the star crosses in the stable region. 



Once the star crosses in the stable region, the amplitudes 
of the three modes damp to zero. The neutron star cools 
via neutrino emission at constant angular velocity until 
the cooling rate becomes slower than the magnetic dipole 
spin-down rate and then spins down emitting magnetic 
dipole radiation. 

The second scenario has Th — Q x 10^ K and is a 
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FIG. 5: Mixed I-II evolution. Th = 6.0 x 10^ K. The 

star cools to Tga = 1.43 in about 10 minutes. At this point 
the dissipative heating of the 3-modes balances the cooling 
and the star spins down, oscillating around its quasi-steady 
solution for ta^b ~ 0.08 yr until it intersects the r-mode sta- 
bility curve at Tgb = 1.25, fib = 0.45. It then cools down at 
constant = f^b until the cooling is again balanced by the 
viscous heating due to the r-mode. The daughter modes do 
not get excited in this part of the evolution. The star then 
spins and cools down in thermal equilibrium on the r-mode 
stability curve for tc-td ~ 1.2 yr until it enters the stable 
region again at Td ~ 0.60, Cld ~ 0.25. 



mixed I-II evolution: the star first reaches thermal evo- 
lution near the high T side of the r-mode stability curve 
(type I evolution) and later, after it intersects the stabil- 
ity curve, the star cools to the low T branch of the curve 
and continues the spin down on this branch (type II evo- 
lution). The trajectory in the — T space is displayed 
in Fig. [H The star cools from Tgi = 10 to Tg « 2 in the 
first minute or so. The r-mode and daughter mode ampli- 
tudes settle in their quasi-stationary states at Tg « 2 and 
subsequently oscillate around these solutions (Fig. [SJa)). 
Increasing shifts the r-mode stability curve to higher 
T and the cooling is stopped at a higher temperature 
than in the T,i = 2 x 10^ K case. At Tga « 1.43 the star 
reaches thermal equilibrium and starts spinning down os- 
cillating around the quasi-steady C = H curve for about 
ta^h « 0.1 yr. Eventually, the oscillations become unsta- 
ble and their amplitude increases until the star intersects 
the r-mode stability curve at (Tgb,fib) = (1.25,0.45). 
By this point the neutron star has lost approximately 
32% of its initial angular velocity, while the temperature 
decreased by about 10% after the star first reached ther- 
mal equilibrium. All points on the high T branch of the 
r-mode stability curve correspond to unstable thermal 
equilibrium (see Appendix D). So, when the trajectory 
reaches this branch of the stability curve, the star cools 
at constant through the stable region until it intersects 
the r-mode stability curve again (assuming the cooling is 
faster than the magnetic spin-down rate). For this evo- 
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FIG. 6: Mixed I-II evolution. Th = 6.0 x 10^ K. Mode am- 
plitudes normalized by their zero-viscosity parametric insta- 
bility threshold shown (a) settling to their quasi-steady states 
in the first minute or so and then oscillate around these solu- 
tions, (b) oscillating in the type I part of the evolution, and 
(c) switching to a type II oscillation as the star cools from the 
low T to the high T branch of the stability curve. 
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FIG. 8: The trajectory of a neutron star in tlie f2 — Tg plane 
is shown for a (1) mixed II-III and (2) type III evolu- 
tion with Th ~ 1.2 X 10^" K. In both cases star cools to 
Tga = 1.82 in about 6 minutes and crosses in the unstable 
region bounded by the low T branch of the r-mode stabil- 
ity curve. In evolution 1 the cooling is balanced by dissipa- 
tive heating close to the low T branch of the stability curve. 
The star settles on the stability curve and spins down on it, 
ta^ci = 6.5 X 10^ sec. As the temperature and spin frequency 
decrease the thermal equilibrium becomes unstable leading to 
growing thermal oscillations. The r-mode passes its paramet- 
ric instability threshold again and the star eventually cools 
(tci^c2 = 1-2 X 10* sec) to the next C = H curve turning into 
a type III evolution. In evolution 2 the star cools at constant 
Cl for ~ 2 hours until Tgbj ~ 1.13. It then spins down for 
th^^c « 1.4 yr on a C = Ji" curve determined by the quasi- 
steady states of all three modes until the trajectory intersects 
the r-mode stability curve again. 



FIG. 7: Mixed I-II evolution. Th = 6.0 x 10^ K. The r-mode 
amplitude for the type II part of the evolution (Tg < Tgpoak) 
is shown as a function of Tg in (a) and as a function of time in 
(b). Part (a) details the initial oscillation of \Ca\ around the 
thermal equilibrium value jC'c,|c=H r-mode. The r-mode even- 
tually settle to its equilibrium value and continues spinning 
down along the r-mode stability curve. Part (b) shows the 
agreement between the r-mode amplitude from the full evolu- 
tion and the thermal equilibrium r-mode amplitude computed 
on the r-mode stability curve |C'a|c=H r-modc(f2cFS, Tg cfs) 
with n between f2c = ~ 0.45 and fid ~ 0.25. 



is now of type II. The settling oscillations of the r-mode 
amplitude are shown in Fig.[7lja). The stability curve acts 
as an attractor. Fig. [7]Jb) plots the r-mode amplitude 
together with its quasi-steady one mode solution as a 
function of time. The agreement is very good until the 
star intersects the stability curve again. At this point the 
r-mode amplitude damps to zero and the star continues 
to spin-down due to the magnetic dipole torque. 

B. Type II 



lution the cooling through the stable region takes about 
tb^c = 14 hours. 

For the r-mode amplitude, the transition between the 
two unstable regions on the different sides of the stability 
curve is shown in Fig. [Bfc). As the star cools through 
the stable region the r-mode amplitude decreases and the 
parametric instability threshold rises. Once the r-mode 
becomes unstable again the star is trapped in thermal 
equilibrium on the r-mode stability curve. The evolution 



The different stages of type II evolutions are mapped 
by the dashed-lined arrows in Fig. [21 In this scenario the 
star reaches the r-mode stability curve, cools through the 
stable region, and settles in thermal equilibrium along 
the low T branch of the r-mode stability curve. These 
evolutions can be approximated by a one-mode model. 

For the Th = 6.0 x 10^ K and h, = 0.60, type II evo- 
lutions occur when Ca(0) < 3.0 x 10"^. They last about 
1.2 yr for B = G. This timescale increases for lower 
B and is independent of the detuning. In these evolutions 
the r-mode amplitude does not grow above its parametric 
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FIG. 9: Mixed II-III evolution. Th = 1.2 x 10^° K. (a) 
The initial r-mode \Ca\ and daughter mode \C^\, |C^| am- 
plitudes. When the evolution reaches thermal equilibrium on 
the r-mode stability curve, the r-mode settles to its one-mode 
equilibrium value |C'a|c=Hono-modo and the daughter mode 
oscillations damp to zero, (b) The r-mode amplitude and its 
thermal equilibrium solution |(7a|c=Honc-modo are shown for 
the type II part the evolution. ICdc^H onc-modc is a good 
approximation for the average r-mode amplitude. 



instability threshold. When Ca(0) > 3.0 x 10"'' mixed 
I-II evolutions occur. These latter evolutions are similar 
to the one in Fig. [5] and follow the same C = H curve. 

We next consider Th = 1.2 x 10^° K. At ^ 0.67, 
initial r-mode amplitudes |Ca|(0) > 5.0 x 10"^ lead to 
mixed evolutions that start as type II and end as type 
III. ICaKO) < 5.0 X 10"^ result in type III evolutions and 
will be described in the next subsection. The trajectories 
for these two evolutions are presented in Fig. [8l 

The transition r-mode amplitude between different 
types of evolutions ICqKO) is very sensitive to both Th 
and r^i. Lowering f2i decreases the driving rate of the r- 



mode. Similarly, increasing Th shifts the r-mode stability 
curve leading to higher viscosity at higher temperatures 
and hence lower growth rate for the r-mode amplitude. 
Lowering the growth rate of the instability makes it more 
likely that the r-mode amplitude will reach its paramet- 
ric instability threshold later. This increases the value 
of the transition amplitude favoring cooling through the 
first C = H curve. 

In mixed II-III evolutions the r-mode amplitude passes 
its first parametric instability threshold close to the low 
T branch of the r-mode stability curve. The dissipative 
heating from the three modes stops the cooling close to 
the low T branch of the r-mode stability curve. Fig. [9l^a) 
plots the initial evolution of the three mode amplitudes. 
As the trajectory of the star in the Q, — T plane converges 
toward the r-mode stability curve, the r-mode amplitude 
settles into its one mode quasi-steady solution and the 
daughter modes damp to zero. Fig. [9l[b) shows the r- 
modc amplitude oscillating around its quasi-stationary 
state. When the equilibrium becomes unstable and the 
oscillations grow, the r-mode passes its parametric insta- 
bility threshold again. The star escapes this equilibrium 
and cools to the next C = H curve in about 4 hours. It 
then spins down on this curve until the trajectory in- 
tersects the r-mode stability curve again. Afterwards 
the mode amplitudes damp to zero leaving the magnetic 
dipole torque to continue to spin down the star. 



C. Type III 

Type III evolutions occur when the star cools through 
the stable region before the r-mode amplitude reaches its 
parametric instability threshold, overshoots the thermal 
equilibrium on the low T branch of the stability curve and 
then settles on the next C = H curve. As an example we 
are examining evolution 2 in Fig. [51 This simulation has 
Th 1.2 X 10^° K,_and the usual f^i = 0.67, Tg; = 10, 
B = 1013 G. Any |C„|(0) < 5.0 x 10"^ for this value of 
Th leads to approximately the same type III evolution. 
The transition amplitude jCaKO) changes with Q.i and Th 
as discussed before. Lower 17^ (higher Th) lead to higher 
transition amplitudes. 

Fig. llOf a) shows the initial r-mode amplitude and 
daughter mode amplitudes settling into their quasi- 
stationary states in the first ^ 40 minutes of the evo- 
lution. The star continues cooling and reaches thermal 
equilibrium at k, 6.5 x 10^ sec. It briefly oscillates 
and then settles on the quasi-steady C — H curve. Fig. 
llOf b) displays the r-mode amplitude \Ca\ and its quasi- 
stationary solution, which coincides with the parametric 
instability threshold amplitude. The agreement between 
the full evolution \Ca\ and its quasi-steady counterpart is 
very good for the first 0.2 yrs or so. As the star cools and 
spins down slowly the viscosity decreases and the star de- 
velops thermal oscillations that are initially unstable with 
growing amplitude. The thermal oscillations become sta- 
ble again when boundary layer viscosity dominates bulk 
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FIG. 11: Runaway evolutions. Th = 6x 10^ K. The r-mode 
amplitude is shown versus time for the two evolutions with 
different viscosity: for the stable case ( ^hb = 2.0 x 10^ sec~^) 
\Ca\ oscillates and settles close to its quasi-steady state while 



in the unstable case (^hb ~ 1.0 x 
grow. Both evolutions have Ca{0) 



10^ 
= 



sec 
1. 



it continues to 



couplings to the r-mode (see Sec lIIIE]) . So, as the model 
becomes more sophisticated and more modes are added 
due to, for example, superfluidity effects, the average de- 
tuning may become several orders of magnitude smaller 
and this scenario could disappear. 
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FIG. 10: Type III evolution. Th = 1.2 x lO" K. (a) The 

amplitudes of the r-mode \Ca\ and of the two daughter modes 
I (7/3 1, \CW are shown settling to their quasi-steady states for 
the first 40 minutes, (b) The r-mode amplitude and its low- 
est parametric instability threshold is shown as a function 
of time for the whole evolution. The quasi-steady solution, 
which coincides with the parametric instability threshold, is 
a good approximation for the average r-mode amplitude with 
an almost exact agreement in the non-oscillatory part of the 
trajectory. 



RUNAWAY EVOLUTIONS 



Mathematically, one can easily see that Fq > Tp + 
T-y corresponds to an unstable evolution by looking at 
the second derivative of the r-mode amplitude evolution 
equation: 



viscosity and the temperature dependence of the viscosity 
changes. For a detuning of 6u!/2n — 10^"* and B = 10^^ 
G, the spin-down lasts about 1.4 yrs. Afterwards, tra- 
jectory of the star intersects the r-mode stability curve 
again. This spin-down timescale is dependent on the size 
of the detuning and the strength of the magnetic field. 
For B < 10^^ G {tb2^e ~ 130) the gravitational radia- 
tion torque dominates the spin-down. When Sut ^ this 
type of evolution no longer exists. The average detuning 
is roughly inversely proportional to the number of direct 



df2 



where Fq = 

jp/Suj, and 



dCg 
df 



— F^ — F-y) 

ra(r;3 + r^)- 



(22) 



4rj 



7Q/(f7|(5J;|) = ■^a/5Lo, ^ 7^/(f2|fc|) 
r7 = 77/(^^1^'^!) = If Fa > Tp 



F^ then the r-mode is unstable no matter how large the 
daughter modes become. 

Assuming a solution of the r-mode amplitude of the 
form Cq oc exp [st] and taking the daughter modes to be 
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FIG. 12: Runaway evolutions. Th = 6 x 10^ K. A similar 
runaway occurs for low cooling. Here /du = 0.001, Ahb = 
5.0 X 10'' sec~^ for both evolutions. The initial amplitudes are 

|C'a|(0),tablo = 2.0 X 10-^ and ICa | (O)unstablo = 2.0 X 10"*. 



In Fig. [TT] we compare two evolutions: one with lower 
viscosity {A\±, = 1.0 x 10^ sec~^) that is unstable and 
the other with slightly higher hyperon bulk viscosity 
(^hb = 2.0 X 10^ sec~^) that is stable. In the unsta- 
ble evolution, + F^ decrease below F^ ~ 'jgh/Soj. 
The amplitude oscillations in the stable evolution damp 
and settle close to the r-mode quasi-stationary solution, 
while, in the unstable evolution, the r-mode amplitude 
continues to grow passing several parametric instability 
thresholds and exciting more inertial modes. 

In our model, for T > Tpcak the hyperon bulk viscosity 
is a decreasing function of temperature. Another way 
to make the r-mode amplitude more likely to overshoot 
its first parametric instability threshold is to lower the 
neutrino cooling. Slower cooling causes the star to spend 
more time at high temperatures and low viscosity. For 
/du — 0.001, evolutions are unstable for jCcKO) > 2 x 
10""' (See Fig. [12] for an example.) As before, when the 
evolution is unstable, F^-|-F-y decreases below F^ and the 
three modes cannot stop the growth of the instability. 



VI. DETECTION OF GRAVITATIONAL WAVES 

We write t he g ravitational wave amplitude ho following 
Watts et al. |4JI as 



5G 



-Er 



(25) 



where d is the distance to the source and the GW fre- 
quency i^GW = 4:1// 3 = (4/3) X r2/(27r). Here v and are 
the spin frequency and angular velocity of the star. 

The gravitational wave amplitude can be written 
in terms of the r-mode amplitude using Eqw = 
-{uJo./m)jGw = 27GiiMi?2|Ca|2fi as 



,2 _ 5GjGRMR^\C^\^n 



(26) 



constant gives the roots 

s± = i(F„-F^-F^)±i {(F„-F^-F^)2 (23) 



+4 



Fa(F^+F^) 



An 



1/2 



- {(F„ - F^ - F^) 



±W(F„+F^ + F^ 



\Cp\^ + \C,\' 



If Fq, < F^ -|- F^ then the s_ mode is always stable and 
the s+ mode is stable only if 



\Cp\^ + \C^\'' <4T^{Tp+T^)n. 



(24) 



For our model 

h, « 1.0 X 10-^^.,^/^, (^^) M^.,Rl,JC^\. (27) 

The r-mode amplitude can be approximated by its 
quasi-stationary solution 

ho « 1.0 X 10--.,^/>f,4i??,5 ('-^) (28) 
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In type I evolutions, we can use the strong viscosity 
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FIG. 13: Type I evolution. Th = 2 x 10^ K. Gravitational 
wave amplitude at d = 100 kpc and spin frequency versus 
time for different values of the magnetic field and of ^hb- 



limit for the quasi-steady mode amplitude 
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7/3 ^ 77 ~ 1/sec is typical at ~ 1 kHz. Note that at 
1/ ~ 1 kHz ^cR ~ 0.05/sec and Suj ~ 1/sec. The gravita- 
tional wave amplitude is independent of detuning and is 
directly proportional to the internal viscous dissipation. 
So, in principle, if fast spinning young neutron stars ex- 



ist and have many modes with high viscosity that are 
near-resonant with the r-mode, their gravitational wave 
signal could provide a measurement of internal neutron 
star physics. 



The type I evolution considered in Sec. IIV Al {T^ = 
2 X 10^ K) is the most optimistic scenario. It has a grav- 
itational wave amplitude of /iq ~ 5 x 10^^^ &i v = 400 
Hz. All values of ft-o reported in this section assume a 
distance to the source of d = 100 kpc unless explicitly 
stated otherwise. The viscous damping rates vary be- 
tween 0.22 sec-i at 900 Hz and 1.26 sec'^ at 300 Hz for 
7/3, and between 0.17 sec~^ and 0.92 sec~^ for 7-y. The 
strong viscosity limit is accurate within 10% of the full 
quasi-stationary solution for v < 500 Hz. This evolution 
is independent of 6lu. 

Fig.[T3]shows the gravitational wave amplitude and the 
spin frequency as a function of time for the Th = 2 x 10^ K 
type I evolution. Lowering the magnetic field can increase 
the duration of the spin-down from 1.1 yr up to about 
13 yr at B « lO^^ G and 14.5 yr at S « 10" G. The 
spin-down timescale is approximately cx A^^ in the limit 
when gravitational radiation spin-down torque is larger 
than the magnetic dipole spin-down. A more detailed 
discussion on how the duration of the instability scales 
with viscosity is presented in Appendix E. 

The mixed I-II evolution with T/^ = 6 x 10^ K (Sec. 
IIV Ap had higher damping rates: 7^ between 2.2 sec~^ 
(at 900 Hz) and 6.7 sec'^ (at 600 Hz) and 7^ is be- 
tween 1.7 sec~^ and 5.1 sec~^ leading to /lo ~ 6 x lO^^'' 
a,t V — 600 Hz. The strong damping limit is accurate 
within 5% for this scenario. The type I part of the evo- 
lution lasts only about a month and in this period the 
spin frequency decreases from 900 Hz to 600 Hz. The 
timescale is independent of the detuning. Although 
is large, the rapid decrease in spin frequency over a rela- 
tively short time would make detection challenging. 

At about V ^ 600 Hz the star intersects the stabil- 
ity curve and cools across the stable region at constant 
spin frequency. It then reaches equilibrium on the low T 
branch of the stability curve. This corresponds to a type 
II behavior and lasts longer (~ 1.2 yr at B = 10^'^ G), 
but leads to a much lower gravitational wave amplitude. 

In type II evolutions the daughter modes are not sig- 
nificant and the spin-down occurs in thermal equilibrium 
Egw — L,j[T), where T — rcFs(^) is determined by 
the CFS stability curve from equating '^gr — "Ja v ■ 

The 

gravitational wave amplitude becomes 



hl = 



5GL,{T) 



GW 



Assuming direct URCA dominates the cooling 
5Gid[/T|/dui?du(T/Te) 



hi 



(31) 
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Taking the square root 



Ti 5GLdufdvRdv{T/Tc 
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3T3 IbGLduUvRMT/Tc) 



(33) 



Plugging in some fiducial values gives 

,-26TiRdv 1 /^lOOkpc 



ho 
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dU 



0.01 I/kHz V d 

1/2 



(34) 
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The gravitational wave amplitude ho in this type of evo- 
lution is independent of detuning and would be a direct 
measurement of the strength of the neutrino cooling. 

In mixed I-II evolutions the type II part of the evo- 
lution lasts the longest, ^1.2 yr for = 6 x 10^ K 
evolution with B = 10^^ G presented in Sec. IIV Al In 
this case Eq. ([M)) gives a gravitational wave amplitude 
decreasing from /iq ~ 5 x 10"^^ when the star first in- 
tersects the r-mode stability curve to /iq ~ 2 x 10"^'' at 
the end of the evolution. The full type II evolution with 
T/i = 6 X 10^ K occurs on approximately the same spin- 
down timescale ~ 1.2 yr and has a similar ho ranging 
between « 6 x 10^^^ and 2 x lO^'^^ for the same distance 
of 100 kpc. 

In mixed II-III evolutions, the type II part of the evo- 
lution is very short. It lasts only about one week for 
the T,, = 1.2 X 10^° K example in Sec. HVB] and leads 
to larger /iq ~ 3 — 5 x lO"^'*. In this short time the 
spin of the star changes by more than 30 % of the initial 
spin frequency, which would make the detection difficult. 
This spin-down timescale for this part of the evolution 
is independent of detuning and also independent of the 
strength of the magnetic field for B < 10^^ G. 

In the limit when 



n\6[d\ 



%v + lp + % - IGR 



» 1 



(35) 



ho « 4.8 X 10-25<,Mi.4i?L.5 ('1^') (36) 



\Cc\ 



\SQ\ 
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The ratio y^jjrh/i'^av + 7/3 + 77 - Igr) is of order 
1. In the type III evolution presented in Sec. IIV CI 
(T/j = 1.2 X 10"^° K ) this ratio varies between 0.5 and 
1.3. The spin-down time changes with detuning, mag- 
netic field and slippage factor. For B = 10^'^ G and 



Sns ~ 0.3, ispin-down 1 yr with /lo ~ 2 X 10 2^ at 
ly — 300 Hz. Lowering the magnetic field to B — 10^^ G 
leads to a spin-down time of ^ 74 yr for the same slip- 
page factor. Lowering the slippage factor lowers spin 
frequency at which the star intersects the r-mode sta- 
bility and can increase the instability timescale signifi- 
cantly. For Sns ~ 0.005 (j// - 80 Hz) and B = lO^^ 
G, the spin-down time is 27 yr. At B = 10^^ G, 
^spin-down ~ 1,400 yr. This timescale also changes with 
detuning. For \6Cj\ « 2 x 10~^ the spin-down occurs on 
a different C = H curve at a lower temperature. For 
~ 0.3 the spin-down time is ^ 75 yr at S = 10^^ G. 
As Slo this scenario disappears. 

An optimistic minimum detectable grav itational wave 
amplitude ho is given by Watts et al. |44| 
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DT, 



obs 



(37) 



where Tobs is the observation time, D is the number of 
detectors and Sn is the power spectral density of the 
detector noise. An integration time of about two weeks 
leads leads to a gravitational wave amplitude of 



2 X 10 



'26 



Sr, 



10-47 

1.21 X 10*^ sec 



1/2 



sec 

1/2 



1/2 



(38) 



^obs 

for Advanced LIGO. The power spectral density Sn of 
Ref. [13 varies between 1.4 x lO""*^ at 300 Hz and 2.2 x 
10"''^ at 900 Hz for their Advanced LIGO broad band 
curve. 

The most optimistic scenarios for gravitational wave 
detection are type I evolutions with ho — 5 x 
10"2'^(lMpc/d)(i//400Hz)2. The gravitational wave am- 
plitude decreases with internal viscous dissipation. For 
a magnetic field oi B = 10^^ G the r-mode instability 
timescale is ^ 1 yr. The timescale increases to ^ 14 
yr when the magnetic field is lowered to B — 10^^ G 
and can be extended further by lowering the dissipation 
(See Fig. ??.) In general, the timescale on which the 
r-mode instability is active is a combination of the grav- 
itational spin-down and the magnetic dipole spin-down 
timescale. The gravitational spin-down timescale due to 
the r-mode is inversely proportional to the viscous dis- 
sipation, while the magnetic dipole spin-down timescale 
is approximately independent of viscosity and oc B~^. 
The An approximate analytic formula for the spin-down 
timescale is derived in Appendix E. 

The core-collapse supernova rate in our galaxy is still 
uncertain. Diehl et al. 50] estimate a galactic core col- 
lapse supernova rate of 1.9 ± 1.1/century from the mea- 
surements of 7-ray radiation from ^^Al. The supernova 
formation rate within a galaxy is typically taken to be 
proportional to the blue band luminosity of the galaxy 
[5lj. Since our galaxy has a blue band luminosity about 
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5 times larger than that of any other galaxy within 100 
kpc, the supernova rate within this distance is approxi- 
mately the same as the galactic rate (See the data from 
Table I of Kopparapu et al. (s^ for blue band luminosity 
of galaxies up to 10 Mpc.) Advanced LIGO may see type 
I evolutions up to 1 Mpc. This distance would include 
about three times as many supernova as in our galaxy be- 
cause 1 Mpc includes the Andromeda galaxy, which has 
more than twice the blue band luminosity of the Milky 
Way [12. 

If a core collapse supernova occurs within 100 kpc while 
Advanced LIGO is on, then the gravitational radiation 
amplitude associated with the r-mode instability ought to 
be large enough to detect in all scenarios we have consid- 
ered. The main complication for data analysis would be 
the relatively rapid, and possibly irregular, neutron star 
spin-down, so that the source is not strictly periodic. In 
the most optimistic scenario we find, we estimate that 
Advanced LIGO would detect gravitational waves from 
r-modes up to a distance of up to 1 Mpc if the neutron 
star was formed less than 15 years ago. No supernova 
has been observed within 100 kpc or even 1 Mpc during 
the past 15 years. Thus, unless a core collapse supernova 
has occurred in a region that is highly obscured from 
Earth, we cannot draw any conclusions from the absence 
of LIGO sources, even in our most optimistic scenarios 
for gravitational radiation associated with the r-mode in- 
stability. Ultimately, if LIGO operates for decades, the 
probability of a core collapse supernova within the galaxy 
is significant. Assuming that our galaxy has a rate of 1 
SN/century, we would see 1 SN/century within 100 kpc 
and 3 SN/century within 1 Mpc. If SN 1987A had pro- 
duced a fast spinning pulsar. Advanced LIGO would only 
be able to see gravitational waves emitted by its n = 3, 
TO = 2 r-mode if the instability timescale is long enough, 
at least 27 years to be detected by Advanced LIGO. For 
our most optimistic scenario, a spin-down timescale of 
~ 30 yr occurs for B = 10^^ G or lower when we decrease 
the hyperon bulk viscosity by a factor of two. However, 
it is unlikely that SN 1987A produced a fast spinning 
neutron star, as electromagnetic energy emitted by the 
pulsar would alter the supernova light curve at late times. 
Observations to date put upper limits on the luminosity 
from a central object in SN 1987A ,53.. .54]. 



VII. LIMITATIONS OF THE MODEL 

There are many ways to improve our calculation. We 
use the eigenmodes and eigenfrequencies of an incom- 
pressible perfect fluid star. These eigenmodes have the 
advantage that they are known analytically and can be 
expressed in terms of Legendre functions labelled by in- 
tegers n and to. This simplifies the computation of bulk"^ 
and boundary layer viscosity and allows us to use the cou- 
pling coefficients computed by Brink et al. [1^ . As stated 
earlier, the inertial modes that comprise the lowest para- 
metric instability threshold change with angular velocity. 



We are using an effective three mode system with typ- 
ical coupling coefficient and detuning. A more realistic 
treatment of the modes would include higher order rota- 
tional effects in the mode frequencies and would model 
the change in frequencies and coupling coefficients with 
angular velocity. This would allow for mode changes in 
the lowest parametric instability threshold and the spin- 
down of the star could be followed more accurately. Ad- 
ditionally, one could include differential rotation [32] and 
mixtures of superfluids or of superfluid and normal fiuid 
regions [ssj . Dissipation rates, particularly from bulk vis- 
cosity, depend on the composition of high density nuclear 
matter which could differ from what we assume. 

We do not include effects due to buoyant forces or mag- 
netic fields in our mode treatment. The magnetic field 
frequency corrections are expected to be small. Morsink 
and Rezania computed this frequency shift Auj/lo 
perturbatively for n = m -I- 1 r-modes by assuming that 
the modes remain unchanged. Following their method 
we compute the shift for several high n inertial modes 
for a constant magnetic field (See Appendix C.) We find 
that for a magnetic field oi B = IQ^^ G and a spin 
frequency of 1 kHz the inertial mode frequency shift is 
Aw/w - B?2(10^^ - 10"^)/i^feff^- This change in fre- 
quency is about three orders of magnitude larger than 
that for the n = 3, to = 2 r-mode, but still much smaller 
than the typical detuning of 5u)/2Q. = lO"'*. So, mag- 
netic field corrections can be ignored. 

Buoyant forces are expected to change the inertial 
mode spectrum at low frequencies. Roughly speaking, 
the Coriolis effect dominates buoyancy when the angular 
velocity of the star 17 is much larger than the Brunt- 
Vaisala angular frequency N . When 2^1 > N the fre- 
quency spacing in the inertial mode sector (inertial mode 
frequencies uj lie between — 217 < w < 20 ) will not be 
affected significantly 35]. Passamonti et aZ. 35] find that 
the Coriolis force is the dominant restoring force for all 
modes when D. > O.^^/nGp. The n = m + 1 r-modes are 
unaffected by buoyancy at all frequencies [s^. As the 
frequency decreases below this value the inertial modes 
behave more like g-modes. As f2 ^ 0, the inertial mode 
frequency is no longer oc O and tends to a constant, lead- 
ing to the loss of resonances with the r-mode. Conse- 
quently jfor a typical Brunt- Vaisala spin frequency of 150 
Hz [3a, |33l the r-mode will lose its near-resonant modes 
roughly when 20 < N {i/ < 75 Hz), and in the absence 
of viscosity it would grow unabated. However, in most 
models the viscous damping rate dominates the gravita- 
tional driving at such low frequencies. For the simple 
viscosity model we consider in this paper the minimum 
spin frequency on the stability curve is higher than 75 
Hz even for low slippage factors ~ 0.01 and a rela- 



ys] Note that incompressible stars do not have bulk viscosity. We 
use an approximation to leading order in the adiabatic index 
that agrees within a factor of two with the n = 1 polytrope 
calculation. See Appendix B for more details. 
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tively high critical temperature for hyperons of Th 10^'^ 
K. See Sec. IIIIE"2l for a more detailed discussion. Thus, 
buoyancy can be safely neglected in our calculations. 



VIII. CONCLUSION 

This paper is the first treatment of the r-mode insta- 
bility that uses a physical model for nonlinear saturation 
in newborn neutron stars. We use one triplet of modes: 
the 71 = 3,771 = 2 r-mode, and two near-resonant inertial 
modes that couple to it. Nonlinear effects become impor- 
tant when the r-mode amplitude grows above its para- 
metric instability threshold. This threshold provides a 
physical cutoff to the r-mode instability by energy trans- 
fer to other inertial modes in the system. 

The behavior we find is richer than expected from pre- 
vious work. We find a variety of evolutions that depend 
on internal neutron star physics, which is still very un- 
certain. The competition between viscous heating and 
neutrino cooling plays the most prominent role. Once 
the cooling is stopped by dissipative heating, the star 
spins down oscillating around quasi-steady thermal equi- 
librium. 

Our model has some schematic aspects. We use typ- 
ical detuning and coupling coefficients that will change 
as the model becomes more sophisticated. The detuning 
is roughly inversely proportional to the number of direct 
couplings to the r-mode, which is a function of the num- 
ber of modes included. Some evolutions are independent 
of this frequency difference between modes. When the 
inertial mode viscosity is larger than the detuning, the 
parametric instability threshold becomes independent of 
the detuning. The threshold in the strong viscosity limit 
coincides with its value for zero detuning among the fre- 
quencies of the mode triplet and still leads to significant 
nonlinear effects. 

Some of our evolutions lead to a gravitational radi- 
ation amplitude that is detectable by Advanced LIGO. 
Assuming a broad band Advanced LIGO curve with two 
detectors and an integration time of about two weeks, 
we find that gravitational radiation would be detectable 
within about 1 Mpc of Earth. The sources would have to 
be young neutron stars within years or perhaps decades 
after formation. The gravitational wave amplitude and 
the duration of the emission depend on the internal dis- 
sipation of the modes of the star. Significant spin-down 
of the star over short periods of time (~ weeks) would 
make detection challenging. However, if fast spinning 
young neutron stars exist and are detected, the gravita- 
tional wave amplitude could give unique information on 
internal neutron star physics. 
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Appendix A: Quasi-Steady Mode Amplitudes 



In terms of amplitudes and phase variables Cj ~ 



|Cj|e"^J Eqs. ^ can be rewritten as 
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where we have defined the relative phase difference as 

The quasi-stationary solution in Eq. (fTT|) is found by 
setting d\Cj\/df — and d(p/df = simultaneously. 



Appendix B: Viscous Dissipation 

This appendix discusses the dissipation we use for the 
r-mode that has mode index j = 4 (ti = 3, tii = 2) and 
two near resonant inertial modes with j = 494 (ti = 14, 
771 — —5) and j — 592 (71 = 15, 771 = 3). Following 
Brink et al. [22l |. we label the hybrid inertial modes so 
that each mode is given a unique number j that is a 
function of principal Legendre index 77, (71 : 2 00), 
azimuthal number ?7i (?7i : — > n — 1) and frequency 
index k (fc : 1 ^ 71 — tti if m 7^ and l^Ti— 771 — lif 
771 = 0). 

^. _ (ri - 1)71(71 + 1) ^ (n - 771 - l)(n - m) ^ ,^ ^ ^^^^ 



We thank Sharon Morsink and Jeandrew Brink for use- 
ful discussions in understanding the effect of high spin 



We compute the bulk viscosity using the modes of the 
incompressible star by taking the dissipation energy to 
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leading order in the adiabatic index Fi 



p 



(41) 



and setting Fi = 2. This ap proximation was proposed 
by Cutler and Lindblom [4^ and adopted by Kokkotas 
and Stergiouluas [49j for the r-mode and by Brink et al. 
EES [2^121 for inertial modes. 

Table H] compares the bulk viscosity timescales for sev- 
eral different inertial modes of n = 1 polytrope computed 
by Lockitch and Friedman (soj with those computed for 
an incompressible model. The difference is small; typi- 
cally about a factor of two or less. For the computations 
in this table we used n-p-e bulk viscosity with a bulk 
viscosity coefficient 



C = 6 X 10 



25 



HzV 

UJ ) 



10 



15 



Tg g cm sec ^ . 



(42) 



TABLE I: Bulk viscosity timescales computed using incom- 
pressible stellar modes (rbuikn=o; used in this work with 
a different bulk viscosity coefficient) and n = 1 polytrope 
(T"buikn=i) inertial modes for an n-p-e gas. The n = 1 poly- 
trope calculation was performed by Lockitch and Friedman 
[301 . Note that, mathematically, the bulk viscosity damping 
rate is zero for incompressible stars. We adopt the attitude 
that the dissipation timescale is computed to leading order in 
Fi and take Fi = 2 (for incompressible stars Fi — * oo). The 
timescales for the two models are roughly within a factor of 
2 of each other. 
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For hyperon bulk viscosity we take the viscous damp- 
ing rates to be 
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The mode integrals /hbj for the inertial modes is com- 
puted using 



1 



Aibj 



F2i?3 



(Fx 



i3 ■ 



V 



(44) 



The r-mode integral is computed by fitting results of Nay- 
yar and Owen T4]. 

We treat the critical temperature Ty^ and the coefficient 
^hb as parameters. Ahb = 5 x lO'^ sec"^ is chosen so that 



the peak of the r-mode stability curve is at roughly 1000 
Hz for all our stable evolutions. We lower it to exhibit 
runaway behavior. 

For the two inertial modes of interest 
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For the r-mode 

-^lib4 = -fhbo 

The relaxation timescale 



0.211. 
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i?hb(r/r^)' 

where we use the reduction factor proposed by Ref. [3! 

^5/4 I ^1/2 

Rhh{TlT^) = exp (^0.5068 - VO.50682 + j/2 

(48) 

w here a = 1 + 0.31182/2, h = 1 + 0.2566^2 and y = 
^1.0-T/r,,(1.456 - 0.157yiyr + 1.764T,,/r). The 
constant is taken to be ti k, 10~^ sec. 

The boundary layer viscosity is computed via Eq. (3) 
of Bildsten and Ushomirsky and is given by 



7, bl(r, 17)= 1,^1(^5/2^2 V 



(49) 



where S-as is the slippage factor and Ij is the mode in- 
tegral for mode j. The slippage factor is the fractional 
difference in velocity between the crust and core of the 
star [13, Sol . Here the constant 



Abi = 3.68 X 10" 



and the mode integral 
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12.5 _i 
sec 
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-^bii" = / dcos 
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(51) 



where is the displacement vector for mode j and 
the density on the boundary layer is taken to be con- 
stant = 1.5 X lO^'' g cm^'^. In these expressions we 
are using the same normalization as Schenk et aZ.(l9j: 
S d^xp = MR^. 

For the r-mode this gives 



7-incomprcssiblc 
^bl4 

Tn=l polytrope 
^bl4 



27r X 10.5, 
27r X 21.8, 



(52) 



where we have used the incompressible value in this pa- 
per. 

For the daughter modes at first two parametric insta- 
bility thresholds the integrals are calculated also using 
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the modes for an incompressible star 

/bi494 ~ 27rx 241.8, /bi 592 » 27r x 140.1. (53) 



Appendix C: Frequency Change due to Magnetic 
Fields 



In this appendix we foUow Morsink and Rezania [SJ] 
to obtain the frequency corrections due to the presence 
of a magnetic field. These corrections are added pertur- 
batively. 

They define a dimensionless magnetic coupling 



KAB 



(54) 



where F is the Lorentz force created by the fluid motion. 
The coefficients can be thought of as the ratio of the 
work done by the perturbed Lorentz force to the total 
magnetic energy stored in the equilibrium star. 

The magnetic coupling coefficients can be written as 

il 



1 



KAB 



4:TTM 



.(1) 
^AB 



J2) 
^AB 



"^B (3) 
^2 '^AB 



(55) 



where 



.(1) _ 

MB ~ 

.(2) _ 

MB — 

,(3) _ 

^AB — 



J d^xS/ X {£_*A^'B)-V X {^B X B) (56) 
J (fx X (V X B) • X (V X B) 
Jd^x (a X B) • (Cb X B). 



We neglect k^^b ^^'^ '^ab- considering a con- 



stant B field and so k 

,(3) 



(2) 
AB 



0. The errors introduced by 

neglecting of the same order as those from ne- 

glecting general relativity [s^ . Assuming the off-diagonal 
entries are small [33 | , the frequency corrections are given 

by 



1 ^ 



(57) 



where M/e is the ratio of magnetic field energy to ro- 
tational kinetic energy. The rotational energy is e = 
MR^rt^ and the magnetic field energy is given by = 

Assuming a constant magnetic field of the form 



Bx ~ -Bo sin a, 

By = 0, 

Bz = Bq cos a, 



(58) 



one can easily compute the magnetic coupling using the 
modes of an incompressible star (Eq. (3.18) together with 



the recursion relations Eq. (A. 1-5) in [24|).') 

For n = m + 1 modes this computation can be per- 
formed analytically [l^ 



_ 47r(TO-f 1)(2to + 3) 

l^r— modes] — X 
+ m — 3 . 2 

1 H sm a 



(59) 



For the n = 3, m = 2 mode 

|kj=4| ~ 22 + 33sin2a (60) 

For the two daughter modes used in this paper 

|k494| ~ 10^(1.8 + 42.0 sin^ a) (61) 
IK592I ~ 10^(2.0 -f 23.5 sin^ a). 

The frequency change is 

M 

t^old ^old 2e 

For the fiducial values used here 



(62) 



_ iQ-12 -^12^12.5 
2e ^ ' 



(63) 



where -812 = 10^^ G. So, the frequency change Aujj/ujj 
for the mode j is 



1^-11^12^12.5 



10" 



(2.1-h3.15sin2a), (64) 



UJ4 ^"1.4''kHz 

^ = 10-9^2:^(1.8 + 42.0sin2a) 

W494 

^ 10-9 -612^^12.5 (2.0 + 23.5 sin^g) 

C^592 -.2 ^ 



The frequency changes are much smaller than the typical 
detuning value of (5w/(20) ~ lO""* for the magnetic field 
values considered in this paper of B < 10^^ G. 



Appendix D: Stability Near Thermal Equilibrium 
One Mode Evolutions 



Consider the one-mode evolution equations 



dt 



2(7Gfl-7a.)|a|^ 



(65) 



.dT 



C{T)^ = 2MR^n\Ca,fja.v-L,{T). 
dt 



We expand each variable to first order around its equi- 
librium value 



T = Te(l-I-CT). 



(66) 
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In equilibrium 



(67) 



This leads to the coupled equations 



dt 

dCr 
dt 



2{jGR-lav)eCa 
^-S'a e'Ta v ^e) 



2T, 



\ dT 



(68) 



TeC{T) 



2E^ 



Ca 



C{T) V dT 



C{T) 



dL^ 



Ct, 



where E^e ^ MR^n\Ca\l- Using Eqs. ^ the coupled 
equations can rewritten as 



dt 
dt 





Ct, 

e 




1 


_ V e 


1 (dL,\ 1 


L^Te) [dTjJ 



(69) 



C(Te) 

We can now write the eigenvalue equation for this system 



A 



51n7Q„ d\nL^\ L^{Te) 



dT dT 



G{Te) 



(70) 



C{Te 



f lav 

\ dT 







with solutions 




2 _ 8£^(re) flav\ 

C{Te) KdTJ, 



Here 



dlnjav 91ni^\ L^{T, 



dT 



dT 



(71) 



(72) 



Points on the right side of the r-mode stability curve T > 
Tpcak have a viscosity with negative slope {d^av / dT)e < 
and are always unstable (one eigenvalue is positive). 
While points on the left side of the r-mode stability curve 
T < Tpcak have {djav/dT)e > are stable if 7e < and 
unstable if 7e > 0. 

In order to gain a better understanding of the thermal 
cycles around the stability curve we write the viscous 
heating U as fraction of the cooling and subsequently 
study the evolution of /. 



U = 2j^,E^^ f{t)L,{T), 



(73) 



where f — 1 corresponds to equal viscous heating and 
cooling and / > 0. We have neglected the viscous heating 
due to the daughter modes as their amplitudes are much 



smaller that that of the r-mode in this scenario. The 
thermal evolution of the system can now be written as 

C(T)^ = U- L,{T) = (/ - 1)L,(T). (74) 

To find the evolution of f we take the time derivative of 
Eq. ([73]). We can then write 



Idf _ 1 dU i f - 1) dL^ 
Jdi ~ U~dt ^ C{T) dT 



(75) 



and 
1 dU 



1 d^^^if -l)L,{T) 



+ 2(7Gfl-7a.)- (76) 



U dt 7„„ dT C{T) 
Plugging Eq. ([7S)l in Eq. (|75|) we can write 

Idf {f-\)L,{T) f 1 dja^v 1 dL, 



f dt C{T) V7c. « dT dT 

+ '^{iGR - lav)- 

Labeling the term in parenthesis as 7 we obtain 
Idf if- l)L,{T) 



(77) 



fdt 



7 + 2(7Gi? I (78) 



where 



7 



C{T) 

1 d-iav 1 dLi, 

jav dT dT 



(79) 



Initially, the star is very hot T ~ 10^° K and cools 
fast: 7 < (the slope of the r-mode stability curve is 
negative for T > Tpcak), / << 1 and jgr > lav So, 
the right hand side of Eq. (|75|) is positive and / grows 
exponentially. All fixed points on the T > Tpcak side of 
the r-mode stability curve are unstable. So, there are 
no evolutions in which the star spins down along the 
right branch (T > Tpcak) of the r-mode stability curve. 
The star can find thermal equilibrium / « 1 in a one 
mode evolution only close to the r-mode stability curve 
IGR ~ lav The thermal oscillations around the stabil- 
ity curve can be understood from this linear perturba- 
tion analysis. Typically, the star overshoots the stability 
curve at first. Once this happens, the right hand side 
of Eq. ([75]) becomes positive (f — 1 and ^gr > lav)- 
This makes / increase and so the star heats. The neu- 
trino cooling (cx T^) eventually balances the heating as 
the temperature increases. This balance happens in the 
stable region igr < lav- At this point (/ = 1 and 
IGR < lav) the right hand side of Eq. ([75)) is nega- 
tive and / starts decreasing (the star cools). The star 
enters the unstable region again and thermal oscillation 
repeats. In this time the angular velocity of the star de- 
creases slowly. So, the next oscillation will have a lower 
amplitude. 

In other words, the thermal equilibrium points at 
fixed fl on the left side of the r-mode stability curve 
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(T < Tpoak; positive slope) are initially stable. The r- 
mode stability curve acts as an attractor. The star ex- 
hibits thermal oscillations at constant angular velocity 
around this curve with the oscillations becoming smaller 
and smaller until the trajectory of the star coincides with 
the r-modc stability curve. As the star spins down the 
viscosity decreases and the heating is slower. If the ther- 
mal equilibrium becomes unstable, then the thermal os- 
cillations restart with growing amplitude until the return 
to thermal equilibrium is no longer possible. The star 
cools until the daughter modes are excited and the vis- 
cosity due to all three modes balances the cooling. Oth- 
erwise, the star continues cooling and spinning down on 
the r-mode stability curve. Eventually, the star enters the 
stable regime again after boundary layer viscosity domi- 
nates bulk viscosity and the slope of the r-mode stability 
curve changes. 



Appendix E: Spin-down Timescales - Scaling with 
Viscosity 



For type I evolutions 



f^^ nc 

''spin — down ^ 



dn-^-r- 



7/377 



(80) 



At first, for simplicity, we set Rhh = Rdu = 1- In this 
case, assuming that hyperon bulk viscosity dominates 



jav oc AhbTg ^jl'', 7^ cx AhbTg ^, 7/3 oc AhbTg ^. (81) 



Using this expression together with Eq. (|84p we obtain 
17, ocA^/^^and t?±_ 



3/2 



« ^hb ana i,rpY„_down « ^hb 



Typically, hyperon bulk viscosity is larger than bound- 
ary layer viscosity for the inertial modes. However, for 
the r-mode boundary layer viscosity can be important at 
low il. For a T/i = 2 X 10^ K evolution without boundary 
layer viscosity with reduction factors included, we obtain 



a scahng of t^f„_down °^ ^hb'^- When boundary layer 
viscosity is included, f2/ becomes fairly independent of 
v4hb changing only by about 1% when Ahb is lowered by 
a factor of two. In this case the spin-down time approx- 
imately scales as 



GR 

spin — down 



«^hb- 



(87) 



In type III evolutions, when we assume that 
1/ tan 02 >> 1, 



spin— down 



On the C — H curve 



(89) 



So, Tg cx Q^/^. On the r-mode stability curve, assuming 
boundary layer viscosity dominates 



Assuming that Tg does not change significantly, we can 
approximate the spin-down timescale by 



^ns V ^ 



(90) 



GR 

spin— down 



1 



OC 



7/377 



"/ -'g/^hb- 



On the C = H curve 



(82) 



(83) 



The spin-down time scales as 

fGR 0-72/41 

spin— down '^ns 



(91) 



When we include the reduction factors in the — 1-2 x 
10^" case we obtain a similar scaling of approximately 

4-GR c-2 
spin — down "^ns ' 



So, 



and 



Tg (X Ali'nf' 



'■spin-down ^ ^hb 



(84) 



(85) 



Assuming the final point is on the r-mode stability curve 
f^^Tg-^Ahb = n"} (86) 



To approximate the total timescale one can use a sim- 
ple interpolation 



ispin-down ~ [l/^spin-down " 



I/^spin— down] i (9^) 



where 



t. 



AID 

spin — down 



3.35 X 10^*5^2 - n--") sec. (93) 



The magnetic dipole spin-down timescale is approxi- 
mately independent of ^hb for type I evolutions. 
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